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Abstract 

This  paper  presents  generic  models  for  the  effect  of  a  chemical 
toxin  on  cells  forming  the  tissue  of  an  organ.  The  models  are 
illustrative,  not  specific  to  organ  or  toxin.  Interactive  response  of 
within-tissue  toxin  and  cells  is  modeled:  cell  capability  to  modify 
(metabolize  or  bind)  toxin  is  represented,  as  is  the  alteration  of  that 
capability  by  toxin  presence.  Both  processes  are  represented  in  the 
context  of  simple  versions  of  the  cell  cycle. 

Toxin-cell  interaction  is  explicitly  represented  in  terms  of  inter- 
cell  chemical  signaling  that  encourages  replacement  or  repair  of 
cells.  The  impact  of  cell  death  by  apoptosis  and  necrosis  is 
represented  in  terms  of  within-tissue  toxin  concentration.  The 
model  is  explicitly  stochastic,  representing  inter-cell  and  toxin  input 
and  within-organ  concentration  in  terms  of  diffusion 
approximations.  Explicit  mathematical  discussion  is  given  of  dose- 
response  function  behavior  at  low  doses. 


1.  Introduction 

Application  of  laboratory  toxicology  data  to  environmental  and  human 
problems  of  risk  assessment  almost  always  requires  extrapolation  of  the  data 
from  the  experimentally-used  dose  regimen  to  the  exposure  conditions  of 
practical  concern,  and  from  the  animal  species  tested  to  the  species  of  concern 
(usually  man).  This  extrapolation  and  estimation  process  is  known  as  chemical 
risk  assessment.  The  risk  assessment  process  has  undergone  considerable 
evolution,  moving  from  a  qualitative  basis  for  decision  making  to  an  increasingly 
quantitative  basis,  and  from  the  use  of  default  assumptions  to  the  application  of 
mechanistic  mathematical  models  as  tools  upon  which  to  base  decisions.  In  the 
context  of  determining  safe  human  exposure  limits  to  potentially  toxic  chemicals, 
there  are  two  sub-tasks  to  be  accomplished:  estimation  of  low-dose  risk  in 
animals,  and  the  subsequent  conversion  of  animal  risk  estimates  to  human  risk 
estimates.  An  authoritative  survey  of  many  of  the  statistical  issues  and 
opportunities  is  given  by  Krewski  and  Franklin  (1991). 

Extrapolation  to  low-dose  effect,  the  first  sub-task,  can  be  accomplished  either 
by  assuming  that  biological  response  varies  in  a  specified  mathematical  manner 
(e.g.  probit,  logit)  with  organ  host  exposure,  or  by  using  physiologically-based 
pharmacokinetic  (PBPK)  models  to  relate  organ  dose  to  host  exposure;  cf. 
Krewski  and  Franklin  (1991),  chapters  8  and  9;  M.  Andersen,  H.  Clewell  Jr.,  and 
C.  Travis  have  been  prominent  in  this  research  area.  In  the  latter  approach,  multi- 
compartment physiological  models  are  formulated  using  actual  tissue  volumes 
from  the  experimental  and  target  species  and  actual  perfusion  rates  to  provide 
for  chemical  transport  between  the  compartments.  Thus  the  pharmacokinetics  of 
high  to  low-dose  extrapolation  become  amenable  to  calculation,  and  external 
measures  of  external  dose  or  exposure  can  be  translated  to  concentrations  in  the 


target  organ  (internal  dose).  Target-organ  chemical  concentration  may  be 
translated  into  estimates  of  risk  if  a  suitable  biologically  acceptable  mathematical 
model  can  be  used  to  relate  chemical  presence  in  the  organ  to  harmful  outcome. 
Such  models,  called  pharmacodynamic,  exist  for  carcinogens,  in  the  form  of  the 
widely-used  linearized  multistage  approximation  for  cancer  dose /response,  and 
multistage  clonal  expansion  models;  cf.  Moolgavkar  (1988).  These  models  are 
based  on  the  general  concept  that  chemical  alteration  of  the  cellular  genes  may 
give  rise  to  permanent,  heritable  changes  in  the  genetic  information  stored  in  the 
cell  nucleus,  and  lead  to  phenotypic  changes  in  the  altered  cells  that  ultimately 
cause  the  formation  of  malignant  tumors.  Simulation  is  sometimes  used  to  study 
stochastic  models  of  carcinogenesis  for  large  numbers  of  cells;  (cf.  Bois  (1992)). 

Unfortunately,  analogous  dose/response  models  are  not  widely  available  for 
toxic  responses  other  than  carcinogenesis.  The  mathematical  expression  of  even 
the  relatively  simple  concept  of  genetic  alteration  leading  to  cancer  involves 
significant  simplification  of  biological  reality,  and  significant  mathematical 
complexity.  This  state  of  affairs  is  exacerbated  when  one  attempts  to  describe  the 
interactions  leading  to  loss  of  cell  function  and  cell  death  in  tissues  of  a  whole 
animal.  The  multi-layered  control  and  response  systems  present  in  an  intact 
living  animal  are  poorly  understood  and  thus  have  not  yet  been  adequately 
modeled.  Nonetheless,  these  control  mechanisms  defend  against  the  majority  of 
toxic  effects  observed  as  a  result  of  chemical  exposure,  and  their  failure 
represents  much  of  what  is  observed  as  the  expression  of  toxicity. 

This  paper  describes  initial  formulations  of  models  for  cellular  response  to 
toxic  insult.  Beginning  with  a  simple  Markov  model,  successive  models  increase 
in  biological  resolution  or  refinement  as  well  as  mathematical  complexity.  Our 
model  formulation  takes  into  account  broad  aspects  of  the  current  state  of 


knowledge  concerning  the  control  and  regulation  of  cellular  properties  in  tissue. 
The  design  of  these  models  is  not  to  provide  a  detailed  description  of  the 
response  of  specific  organs  to  toxic  insult,  but  rather  to  begin  the  process  of 
mathematically  describing  recent  significant  advances  in  knowledge  that  have 
been  made  in  cell  biology,  hormonal  regulation,  and  hormetic  control  mech- 
anisms that  regulate  cell  response  to  toxicants.  Some  modifications  of  the  model 
formulations  may  be  needed  to  apply  them  to  specific  organs;  for  one  thing  the 
complex  geometry  of  organ  architecture  is  not  represented.  To  have  included 
complex,  three-dimensional  relationships  between  different  cell  types  would 
have  increased  the  mathematical  complexity  considerably.  But  even  in  its  present 
rudimentary  state  the  model  is  useful  for  stimulating  discussion  and  suggesting 
experimentation,  with  the  goal  of  furthering  the  understanding  of  cellular 
response  to  toxic  substances  including  the  effects  of  tissue-mediated  responses 
(e.g.  mitogenic  stimulation).  The  current  formulation  of  the  model,  or  extensions 
thereof,  is  suitable  for  experimental  evaluation  in  cell  culture  if  the  experimental 
conditions  are  properly  chosen  and  chemical-specific  details  are  added. 

We  provide  the  following  abbreviated  biological  background.  The  tissue 
making  up  organs  almost  always  involves  a  complex  geometrical  juxtaposition  of 
several  cell  types  having  specific  and  often  overlapping  functions.  This  tissue 
architecture  may  be  maintained  by  the  presence  of  a  nonliving  matrix  of  proteins 
collectively  known  as  a  basement  membrane.  Interactions  between  this  membrane 
and  the  cells  are  known  to  be  critical  to  its  stability  as  a  mature,  functioning 
entity.  The  whole  of  the  tissue  is  permeated  by  branching  blood  vessels,  each 
generation  of  which  is  successively  smaller;  these  serve  to  provide  a  constant, 
stable  milieu  in  which  the  cells  exist.  Nutrients,  control  signals  in  the  form  of 
specific  biomolecules,  and  xenobiotic  chemicals  are  brought  to  the  immediate 


surroundings  of  the  cell  by  the  vasculature  in  the  tissue  and  cellular  metabolic 
products;  the  products  of  energy  metabolism  and  cellularly  derived  control 
biomolecules  are  removed  from  the  cellular  microenvironment  by  the  same 
means. 

The  state  of  the  cell  at  any  time  is  a  reflection  of  its  age,  the  summation  of  the 
control  chemicals  reaching  and  leaving  the  cell,  the  effects  of  xenobiotic 
chemicals  present,  if  any,  and  the  state  of  the  cells  surrounding  it.  Cellular 
contact  with  the  surrounding  cells  and  basement  membrane  act  to  supply 
chemical  signals  that  modulate  its  activity.  A  given  cell  may  be  (i)  nearly 
quiescent,  (ii)  active  biochemically,  i.e.,  producing  metabolites  of  absorbed 
materials  for  its  own  use  or  for  export,  (iii)  in  a  state  of  stress  due  to  shortage  or 
excess  of  biochemical  molecules,  (iv)  in  a  process  of  programmed  cell  death 
(apoptosis),  (v)  in  the  process  of  dying  from  chemical  insult  (necrosis),  or  (vi) 
dividing  to  form  progenitor  cells  in  response  to  a  need  to  replace  cells  already 
lost;  these  conditions  need  not  be  mutually  exclusive  at  any  point  in  time.  Some 
specific  cells  may  alter  or  completely  change  their  observable  characteristics  in 
response  to  chemical  signals.  The  most  well-known  of  these  cells  are  the 
pluripotent  stem  cells  of  the  hematopoietic  system. 

The  models  to  be  presented  here  represent  some,  but  by  no  means  all,  of  the 
features  mentioned.  They  provide  a  basis  for  proceeding  further. 

2.  Model  Structure 

In  the  following  sections  of  this  paper  we  present  a  sequence  of  mathematical 
models  that  represent  interactions  between  the  cells  in  an  organ's  tissue  and  a 
toxic  chemical  or  agent  (hereafter  called  toxin)  originating  externally  and  entering 
the  organ  and  coming  into  contact  with  the  cells,  possibly  in  modified  form.  The 
initial  models  postulate  a  mature  organ  with  capacity  to  contain,  at  most,  a  fixed 


number  of  cells.  Broadly  speaking,  the  organ's  cells  attempt  to  remove  (e.g. 
metabolize  or  bind)  the  toxin,  but  are,  to  a  dose-dependent  degree,  also  affected 
by  the  toxin  concentration.  We  distinguish  between  cells  that  are  in  one  of  two 
stages  of  the  cell  cycle:  these  we  call,  oversimplifying  the  true  cell  cycle, 
functional,  and  dividing  (S-phase  and  mitotic).  Functional  cells  are  assumed  to  be 
capable  of  removing  toxin,  but  are  also  susceptible  to  premature  death  because  of 
the  toxin's  action.  Dividing  cells  are  in  the  process  of  DNA  replication  and  actual 
mitosis,  here  splitting  into  two  daughter  cells;  while  in  this  stage  the  cells  are 
assumed  incapable  of  removing  toxin.  However,  toxin  presence  is  assumed  to 
affect  the  cell  cycle,  e.g.  by  shortening  the  functional,  and  lengthening  the 
dividing,  stages.  The  model  proposed  emphasizes  inter-cellular  signaling:  when 
a  cell  dies,  i.e.  ceases  to  function,  it  effectively  requests  one  other  functioning  cell 
to  divide  to  produce  its  replacement.  It  is  biologically  plausible  that  signaling 
occurs  only  to  neighboring  cells,  but  the  present  model,  being  without  spatial 
characteristics,  does  not  literally  recognize  this  restriction.  Also,  signals  are 
assumed  to  be  specific:  one  signal  is  emitted  per  cell  death,  so  eventually  one 
replacement  occurs.  Alternative  formulations  are  sometimes  plausible  and  are 
easily  studied. 

3.  Single-Stage  Markov  Model 

We  now  present  several  alternative  simple  Markov  models  for  cell-toxin 
interaction,  emphasizing  the  representation  of  signaling.  Although  cells  can  age 
and  exist  in  different  stages  of  functionality  we  omit  an  account  of  this  at  first, 
but  return  to  it  later. 

Let  D(t)  represent  the  number  of  functional  cells  in  the  organ  at  time  t,  and 
M(0  the  number  of  dividing  (mitotic)  cells.  If  a  functioning  cell  that  dies 
immediately  induces  another  to  begin  division  (instant  signaling)  then 


Co  =  D(t)  +  2M(t)  (3.1) 

where  Co  is  the  maximum  organ  size,  here  assumed  constant;  realistic 
elaborations  are  possible.  Note  that  organ  size  is,  strictly  speaking,  D(t)  +  M (f), 
which  is  the  number  of  cells  present  at  time  t;  Co  acknowledges  that  space  for  all 
cells  must  leave  room  for  each  cell  currently  in  mitosis  to  eventually  become  two 
differentiated  cells.  We  describe  {D(0,  M(f)}  as  a  birth  and  death  process  in 
continuous  time,  with  D(t)e  (0, 1,  2, ...,  Co).  Typically  Co  is  large,  i.e.  of  order 
1010-1012  for  mature  humans,  and  this  will  be  exploited  to  carry  out  an 
asymptotic  analysis. 

Let  T(0  denote  the  concentration  of  toxin  in  the  organ  at  time  t.  We  represent 
{T(0}  as  a  diffusion  process  whose  drift  and  diffusion  coefficients  are  both 
influenced  by  D(t)  so  as  to  represent  metabolic  or  binding  action,  hence  toxin 
removal,  and  by  the  toxin  input  to  the  organ,  xtf),  as  well  as  the  instantaneous 
concentration,  T(t),  itself. 

In  the  following  we  describe  specific  models,  beginning  with  the  most  simple 
and  (potentially)  least  realistic,  but  also  the  most  parsimonious.  We  explicitly 
include  the  effect  of  inter-cellular  signaling. 

3.1  Cellular  Signaling 

Suppose  a  cell  dies  from  some  cause:  apoptosis  or  necrosis.  The  effect  of  its 
signaling  to  other  cells  for  replacement  can  be  captured  by  postulating  that  a 
ghost  (release  of  growth  factor)  of  the  newly-dead  cell  now  exists,  the  function  of 
which  is  to  find  and  interact  with  a  neighboring  live  cell,  inducing  that  cell  to 
divide.  The  ghost  then  disappears  after  its  purpose  is  served.  One  can  thus 
represent  the  delay  in  mitotic  response  to  cell  death  by  adjusting  the  ghost's 
search  rate  for  functional  cells. 


6 


The  following  describes  a  possible  Markov  generator  that  includes  the 
signaling  (ghost)  influence.  Let  G(t)  denote  the  number  of  signals  present  in  the 
organism  at  time  t,  and  introduce  the  constraint 

Co  =  D(r)  +  G(0  +  2M(0.  (3.2) 

The  following  then  describes  the  possible  transitions  allowed  to  occur  in  a  small 
time  interval,  (t,  t  +  dt). 

Signaling  Model 
(Markov  Generator) 


t 
D(t),G(t),M{t); 

IX*) 


t  +  dt  Probability 

-4  D(t)  - 1,  Git)  +  1M0  X*D(t)dt 

T(t)  +  8*dt  +  odW(t) 
(cell  death,  signal  creation) 

->  D(f)  - 1,  Git)  -  1,  Mit)  +  1  6  *G(t)(D(t)/C0)dt 

7X0+  6*dt  +  odW(t)  (3.3) 

(cell  signaled  to  divide,  signal 
disappears) 

-»  D(t)  +  2,  G(0,  M(f )  - 1  n*M(t)  dt 

T(t)  +  8*dt  +  odW(t) 
(two  newly-divided  cells  emerge) 


All  rates  of  cell  transition,  namely  A*,  & \  and  (i  are  presumed  to  be  influenced 
by  current  toxin  concentration,  7^0-  In  turn,  the  latter  is  influenced  by  the 
metabolic  or  binding  capability  of  the  cells  to  remove  toxin;  the  net  mean 
increase  of  toxin  in  the  short  run  is  6*(T(t),D(t),x(t)).  Our  subsequent  analysis 
does  not,  in  principle,  depend  upon  specific  functional  forms  for  any  of  the  above 
parameters;  adequate  smoothness  and  differentiability  is  assumed  where 
necessary.  When  specific  functions  are  required  we  shall  use  these: 


£(T(t))  =  x/T{t)  (3.4,a) 

n\T(t))  =  n0e-»T{i)  (3.4,b) 

e*(T(t))=60e-e'T{t);    *>V>0.  (3.4,c) 

Note  that  it  has  been  assumed  that  increased  toxin  level  increases  cell  death  rate, 
decreases  the  rate  of  completing  mitosis,  and  slows  the  rate  at  which  ambient 
signals  reach  their  destinations.  The  net  effect  is  to  reduce  the  number  of 
functional  cells.  All  of  the  above  dependencies  are  hypothetical  and  illustrative 
only.  In  particular 

T(t)(D(t)/C0) 


8*(T{tW)At))=  *M-v-^{5 


no 


(3.5) 


where  in  the  latter  T*(0  represents  toxin  input  rate  to  the  organ,  and  the 
remaining  term  v  (T(f)D(f)/Co)/  1  +  k*  T(t)\  represents  the  toxin  concentration 

rate-limited  (Michaelis-Menton)  reduction  by  functional  cells. 

It  has  been  observed  that  the  realistic  effect  of  toxin  presence  on  A*,  /i*,  0*,  etc. 

at  time  t  may  involve  the  entire  past  history  of  such  exposure,  e.g. 

A  (t)  =  X  (T(x),a<x<t);  the  explicit  function  may  be  illustrated  by  the  form 


A*(0  =  exp 


j'/(T(x)dx) 


A  further  model  enhancement  would  recognize  that 
toxin  presence  in  the  organ  may  result  in  cells  that  complete  completing  mitosis 
giving  rise  to  damaged  cells,  e.g.  possessing  DNA  adducts.  These  cells  now 
become  susceptible  to  repair. 
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3.2  Differential  Equations  For  Deterministic  Approximations 

Let  Co  denote  the  number  of  spaces/holes  that  cells  may  occupy.  Then 

Co  =  D(f)  +  Git)  +  2M(t)  (3.6) 

since  the  cells  in  division  potentially  require  two  holes.  Assume  that  if  Co  »  1 
then 

D(f)/C0  -»  ait),  Git) /Co  -*  yit),  IXfl/Co  -*  Pit)  (3.7) 

in  probability  as  Co— »°°  where  the  latter  functions  are  O(l).  The  following 
differential  equations  result  from  direct  manipulation  of  (3.3): 

£21  =  -Xa{t )  -  6a(t  MO  +  2/i(l  -  a(()  -  y(())| ;  (3.8) 

the  last  term  occurs  by  virtue  of  (3.6).  Also, 

^P  =  Xa(t)-6a(t)y(t)  (3.9) 

at 

*l  =  5  =  T-^-  (3.10) 

dt  1  +  Kp 

We  re-emphasize  that  all  parameters  are  implicit  functions  of  (current)  toxin 

level,  T{t),  and  possibly  also  absolute  time,  t,  which  can  be  viewed  as  the  age  of  a 

mature  organ,  rather  than  an  individual  cell.  This  latter  step  is  not  taken  in  this 

paper. 

Note  that  in  (3.8),  (3.9),  (3.10)  the  original  parameters  of  (3.3)  must  be  scaled: 
KiTit))  is  replaced  by  Xifiit)),  e*(7X0)  by  0(#O),  /i*(T(r))  by  nipit)),  v\Tit))  by 
viPit)),  t*  by  tCo,  and  k  by  kCq. 

In  Section  6  such  scaling  is  exploited  more  extensively  to  deduce  stochastic 
behavior  of  system  state  values. 


3.3  Solution  Without  Toxin  Input 

Suppose  no  toxin  input  exists,  so  T(t)  =  0,  Vf  >  0.  Then  a  steady-state  solution 

to  (3.8)  and  (3.9)  may  occur:  set  derivatives  =  0  and  solve  the  resulting  equations 

to  find 

y  =  \,       provided    0<-<l  (3.11) 

6  6 

«  =  4f^}.  (3,2) 

The  latter  implies  (i)  that  the  larger  the  signal  activity  rate,  6,  the  smaller  the 
ambient  signal  population,  and  the  more  quickly  does  signaling  induce  another 
cell  to  begin  division;  furthermore  (ii)  if  8  »  A,  the  death  rate,  signifying  highly 
efficient  signaling,  then  a  — » n/(2k  +  /i),  which  is  equivalent  to  thinking  of  cells  as 
behaving  in  pairs:  once  division  completes  the  two  daughters  effectively  compete 
to  die  (rate  2A). 

Although  the  equations  (3.8)  and  (3.9)  are  non-linear  and  apparently  cannot 
be  solved  explicitly  one  can  obtain  an  approximate  time-dependent  solution  as 
follows: 

(a)     assume  ){t)  quickly  reaches  steady  state,  so 

y(t)  =  X/e;     0<t  (3.13) 

This  is  a  standard  and  useful  approximation  technique  often  invoked  in 
biochemical  kinetics  problems  that  is  variously  abbreviated  the  quasi-steady-state 
assumption  (QSSA)  or  the  pseudo-steady-state  hypothesis  (PSSH);  see  Segel  and 
Slemrod  (1989)  for  a  careful  exposition  of  its  rationale.  The  QSSA  approach  gives 
rise  to  the  classical  Michaelis-Menton  used  widely  in  pharmaceutical  kinetics  and 
compartment  models;  it  already  appears  implicitly  in  our  (3.5)  and  (3.10). 
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(b)    introduce  (3.13)  into  (3.8)  to  find  the  approximation  a(t ) : 

da{t)/dt  =  -(2A  +  n)a  +  n(l  -  X/6), 
the  solution  of  which  is 

&(t)  =  a(0)e~{n+^  +  [m(1  "  ^)/(2A  +  /i)](l  -  e<ZX+^ ) 
as  in  (3.12). 


(3.14) 


(3.15) 


3.4  Dose-Response  for  Small  Steady  Toxin  Exposure 

Of  interest  in  risk  analysis  is  the  behavior  of  the  dose  response  curve  for  small 
values  of  (toxic)  dose.  We  approximate  this  by  finding  expressions  for 


da(r) 


dr 


dy(r) 


t  =  0'    dr 


and-^ 
t  =  0  dr 


T  =  0" 


Suppose  tit)  =  t  =  t  ~  0,  a  constant,  and  that  exposure  has  proceeded  for  some 
time  so  that  a  steady-state  condition  has  been  reached;  this  is  modeled  by  letting 

M0=M0  =  ^  =  o  iR  (3g)_  (310)   To  stud     the  dose-response  for  r«0, 
dt         dt       dt  j  r 

differentiate  (3.10)  with  respect  to  t  using  the  function  S  in  (3.5)  (drop  subscript 

"o"  for  convenience) 


dvdp    ap       vda     p  va    dp  |      vccp     dp 

dpdrl  +  icp     dr  1  +  kP     1  +  tcp  dr     (1  +  tcp)2  dr 


(3.16) 


Since  p  =  0  for  t  =  0,  equation  (3.6)  yields  for  the  (scaled)  rate  of  toxin 
concentration  at  low  exposure 


*■* 


t  =  0     va 


r  =  0     v^(i-|)  0 


(3.17) 


This  shows  explicitly  how  low-dose  toxin  level  in  the  organ  tissue  increases  as  a 
function  of  cell-cycle  parameters  X  and  \i,  signaling  efficiency,  6,  and  toxin 
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metabolic  or  binding  rate  v.  All  of  these  parameters  are  evaluated  at  very  low 
(zero)  toxin  levels. 

Next  consider  steady-state  rate  of  change  of  scaled  functional  cell  fraction, 

da 

— -  =  ax  =  a  a  •  pr;  this  is  the  slope  of  a  dose-(functioning  cell  count)  response 

ox  H 

curve  at  very  low  levels  of  toxin.  If  (3.12)  is  differentiated  and  some 
rearrangements  made,  and  if  y=  XI 6  as  in  (3.11),  then 

H  =hV(2  +  A/M  +  2(l-y))  +  ^(2A(l-y))  +  ^(2A  +  /z)(/z/g)y]       ^     2  .   (3.18) 

Since  /?T  >  0,  the  sign  of  Ox  is  determined  by  the  above  linear  combination  of  the 
derivatives  of  the  cell  cycle  and  signaling  parameters  A,  fi,  and  6  at  x,  and  hence 
(3,  equal  to  zero.  Note  that  if  the  illustrative  relations  (3.4)  are  invoked  as  given, 
then  ap  <  0,  hence  a?  <  0  so  a  small  increase  in  toxin  induces  a  decline  in  the 
number  of  functioning  cells;  this  is  intuitively  acceptable.  However,  a  change  in 
sign  of  any  of  the  parameters'  derivatives  is  capable  of  reversing  the  sign  of  the 
slope  of  the  dose  response  curve,  thus  representing  a  tendency  for  hormesis  at 
low  dose. 

4.  Models  with  Cell  Aging 

In  this  section  we  generalize  the  models  of  the  last  section  to  allow  cells  to 
age,  and  hence  exist  in  possibly  different  stages  of  functional  effectiveness. 

We  restrict  attention  to  just  two  age  stages,  which  we  term  new  (n)  and  old  (o), 
distinguishing  them  by  subscript.  An  arbitrary  number  of  such  age  states  can  be 
introduced,  but  at  the  cost  of  increasing  the  number  of  parameters. 

Let  Dn(t)  (respectively  D0(f))  represent  the  number  of  new  (respectively  old) 
cells  at  time  t.  M(t)  represents  the  number  of  dividing  (mitotic)  cells  at  time  t  and 
G(t)  the  number  of  ambient  signals  or  ghosts  as  before.  The  spatial  constraint  is 
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Co  =  Dn(t)  +  D0(t)  +  Git)  +  2M(f). 
Here  are  the  transition  probabilities  assumed: 

Signaling  Model  with  Age  Dependence 
(Markov  Generator) 


(4.1) 


t  t  +  dt 

Dn(t),Doit),    ->         Dn(t)-l,D0(t)tG(t)  +  lM(t) 

Git),  Mit);  Tit)  +  8*dt  +  odWit) 

T(t)  (new  cell  death,  signal  creation) 

->         Dnit),D0it)-l,Git)  +  l,Mit) 
Tit)  +  8*dt  +  odWit) 
(old  cell  death,  signal  creation)l 

->      Dnit)  -  1,  DM  Git)  -  l,M«)  +  1 
Tit)  +  8*dt  +  odWit) 
(new  cell  signaled  to  divide, 
signal  disappears) 

->     Dnit)f  D0it)  - 1,  G(t)  -  1,  M(0  +  1 
Tit)  +  8*dt  +  odWit) 
(old  cell  signaled  to  divide, 
signal  disappears) 

Dnit)  +  2,D0it),Git),Mit)-l 

Tit)  +  8*dt  +  odWit) 

(two  newly-divided  new  cells 

emerge) 

->        Dnit) -h  Doit)  +  I,  Git),  Mit) 
Tit)  +  8*dt  +  odWit) 
(a  new  cell  becomes  an  old  cell) 


Probability 
£nDnit)dt 


X0D0it)dt 


dnGit)iDnit)/C0)dt 


eoGit)iD0it)/Co)dt      (4.2) 


tfMit)  dt 


<t>*Dnit)dt 


Once  again  all  rates  of  cell  transition,  namely  aJ,,  X0,  6n,  60,  /z*  are  presumed 
to  be  influenced  by  current  toxin  concentration.  In  turn,  the  latter  is  influenced  by 
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the  metabolic  capability  of  the  cells  to  remove  toxin,  5*(T(0,  Dn(t),  D0(t),  tit)). 
When  specific  functions  are  required,  we  shall  use  for  illustration 

5-(m(wwo)=Afv;M-v;«  (4.3a) 

±  +  KnT(t)  1  +  K0l{t) 


^(r(0)=Vo 


re^T®-l 


+1  > 

J 


J 


p+q  =  l,  p>0 


4(r(0)=^,o^',T(0 


(4.3b) 

(4.3c) 
(4.3d) 


ti*(T{t))  =  tioe-^T^  (4.3e) 

The  form  of  (4.3b)  permits  the  adjustment  of  organ  toxin  concentration-caused 
cell  death  rate  to  be  flexibly  adjusted;  if  Xn\  <  0  the  inner  bracket  must  have  its 

sign  reversed:  (    )+.  As  before,  the  above  forms  are  hypothetical. 

4.1  Differential  Equations  for  Deterministic  Approximations 

Let  Co  denote  the  number  of  spaces /holes  that  cells  may  occupy.  Then 


Co  =  Dn(O  +  Do(O  +  G(0  +  2M(0. 


Assume  that 


Co  <-0  M)  L0 


(4.4) 


in  probability  as  Co  — >  °°-  The  following  differential  equations  result  from  direct 
manipulation  of  (4.2)  with  scaled  rate  functions 

^^  =  -[AB+0]an(O-^a„(Oy(O  +  4l-«n(O-ao(O-y(O]  (4.5) 

fSM  =  -X0a0(t)  +  0an{t)  -  eoa0(t)y(t)  (4.6) 

at 
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Ml = -[e„a„(() + e0a0(t)]y(t) + W) + V*„(<)  <47> 

These  equations  can  be  solved  numerically,  but  not  in  closed  parametric  form. 

4.2  Solution  Without  Toxin  Input 

Suppose  no  toxin  input  exists  so  T(t)  =  0,  t  £  0.  Then  a  steady-state  solution  to 
(4.5)  -  (4.7)  may  occur;  set  derivatives  equal  to  0  obtaining  the  following 
equations. 

0  =  -(A„  +  <f>)an  -  6nany  +  /i[l  -an-a0-y]  (4.9) 

0  =  -Xoao  +  <t>an-6oaoY  (4.10) 

0  =  [-Bnan  -  eoa0]y  +  ^  +  X0a0.  (4.11) 

Equation  (4.10)  yields 

a„=(l/0)[Ao  +  0oy]ao.  (4.12) 

Substituting  (4.12)  into  equation  (4.11)  yields 

o  =  [-dn(y<p)[x0 + e0y]  -  e0]r + {K[{V<t>\K + Ml + A*)  (413) 

after  dividing  through  by  ocq. 
Rewriting  (4.13)  we  obtain 

o = (-ene0/<i>)Y2  -  y[en(y<t>)x0 +e0-  K{Vt%] + [( V#) + i]V      (4.14) 

Thus  y  satisfies  a  quadratic  equation.  Let  y0  be  the  positive  solution  to  the 
equation.  If  An(0)  =  0,  as  in  (4.3b) 


{k^Hk^)2+*ki 


yo=_U_* nj  ,V     w ^ ^  (415) 
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If  y0  <  1,  then  y=  y0.  If  y0  >  1  then  ambient  signals  and  not  functional  cells 
dominate  the  organ,  which  in  reality  would  be  long  defunct  even  without  toxin. 

If  —r  <  1  and 
0 


er 


.     0. 


.      0J 


(4.16) 


then  y0  <  1. 

To  find  olq  equation  (4.9)  is  used.  Putci  =  (l/#)Uo  +  Oofi.  Then 


/41  ~  To]  =  {%  +  0)q  +  0»tt  +  /4ci  +  l]}a0  =  {2A0  +  //[q  + 1]  +  2A„c1}or0. 


Thus, 


and 


«o  =  M[l "  y0]/{2A0  +  /i[l  +  q]  +  2  Vl} 


a„=cia0. 


(4.17) 


(4.18) 


4.3  Dose-Response  For  Small  Steady  Toxin  Input 

As  mentioned  earlier,  the  behavior  of  an  organ-level  dose-response  curve  for 
small  values  of  toxic  dose  is  of  interest  in  risk  analysis.  In  this  section  we  indicate 
that  the  derivative  of  various  cellular-level  responses  with  respect  to  t  may  be 
evaluated  at  x  =  0. 

If  t  =  0  then  P  =  0.  Set  the  left-hand  side  of  (4.8)  equal  to  zero,  then 
differentiate  the  right  side  with  respect  to  x.  After  setting  p  -  0,  this  yields 


0  =  l-[vnan  +  voao]-£ 


Thus 


dp 
dx 


t  =  0     vnan+v0a0 


>0, 


(4.19) 


(4.20) 


as  is  physically  plausible. 
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It  may  be  shown  by  further  differentiation  that  under  the  specific 
parameterization  of  (3.4),  or  suitable  alternatives,  all  derivatives  of  long-run 
cellular-level  response  can  be  explicitly  evaluated  at  t  =  0,  as  was  shown  earlier. 
The  complex  details  are  omitted  here. 

5.  A  Model  with  Cell  Aging  and  Organ  Tissue  Growth 

In  this  section  we  generalize  the  previous  model  to  allow  cells  to 
spontaneously  enter  mitosis  and  multiply,  thus  effectively  increasing  the  organ 
size,  where  the  latter  is  defined  as  the  number  of  cells  in  existence,  in  either  state. 
In  the  previous  section  cell  replication  only  occurred  to  replace  other  cells  that 
previously  have  died,  possibly  as  a  result  of  toxin  action. 

Consider  the  region  to  be  occupied  by  cells  to  consist  of  spaces /holes  that 
may  be  filled  by  cells.  Let  Co  be  the  maximum  number  of  such  holes,  in  this  case 
the  number  in  an  essentially  mature  organ.  Initially  not  all  such  holes  are  active, 
for  the  organ  is  immature  and  hence  growing.  Hence  let  C(t)  be  the  number  of 
active  holes  at  time  t ,  while  Co  -  C(t)  are  currently  quiescent  or  inactive.  Now, 
generalizing  the  previous  setup,  let  the  current  number  of  active  holes  /spaces  be 

C(t)  =  Dn(t)  +  D0(t)  +  Git)  +  2M(t)  (5.1) 

and  augment  the  previous  model  to  allow  C(t)  to  gradually  grow,  although 
remaining  bounded  by  Co- 

5.1  Organ  Growth  Stimulated  by  Signaling 

In  this  model  the  presence  of  inactive  holes /spaces  effectively  encourages  the 
active  cells,  Dn(t)  and  D0(t),  to  increase  the  number  of  signals  and  hence  the 
number  of  active  spaces  capable  of  accomodating  cells.  We  are  implicitly 
modeling  the  competing  effects  of  positive  and  negative,  or  inhibitory,  growth 
factors. 
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The  model  is  specified  by  the  following  transitions. 


Signaling-Driven  Organ  Growth 
(Markov  Generator) 


t  t+dt 

Dn(t\  D0(t),    -+     Dnit),  D0it),  Git)  +  1,  Mit),  Cit)  +  1 
Git),  Mit),  Tit)  +  Sdt  +  odWit) 

dt)  Tit)  (a  quiescent  cell  space  becomes 

active,  signal  creation) 

-»     Dnit)  - 1,  D0it),  Git)  +  1,  Mit),  Cit) 
Tit)  +  5dt  +  odWit) 
(new  cell  death,  signal  creation) 

->      Dnit),  D0it)  - 1,  Git)  +  1,  M(0,  C(0 
Tit)  +  6dt  +  odWit) 
(old  cell  death,  signal  creation) 

->  D„(f)  +  2,  Do(0,  G(0  +  1,  M(0  -  1,  C(0 
Tit)  +  6dt  +  odWit) 
(two  newly-divided  cells  emerge) 

-»     D„(f)  - 1,  Do(0  +  1,  G(0,  Mit),  Cit) 
Tit)  +  8dt  +  odWit) 
(new  cell  becomes  old) 

->  Dnit)  -  1,  D0(f),  Gft)  - 1,  M(0  +  1,  Cit) 
Tit)  +  Sdt  +  odWit) 
(new  cell  signaled  to  divide, 
signal  vanishes) 

->  Dnit),  D0it)  -  1,  Git)  - 1,  M (0  +  1,  C(0 
7X0  +  Sdt  +  orfW(0 
(old  cell  signaled  to  divide, 
signal  vanishes) 


Probability 


[&>b(0+£do(0] 


i- 


c(0 


A 


^D„(0<if 


l0D0it)dt 


H*Mit)dt 


<fDnit)dt 


0nG(t) 


C(t) 


dt 


e0G(t)^-dt 

0    K  }  C(t) 


(5.2) 
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Once  again,  all  rates  of  cell  transition,  namely,  A„,  A0,  €t,  0*,  \i,  £„,  <f;0  are 
presumed  to  be  influenced  by  current  toxin  concentration.  When  specific 
functions  are  required  we  shall  use 


&(T(t))=S°M<<im 


in  addition  to  those  of  (4.3  a-e). 


(53) 


Once  again  assume  that  as  Co  -*  <*° 

L0  Lo  c0 

in  probability. 


m 


^-*m^-*m  (5.4) 


The  following  differential  equations  result  from  direct  manipulation  of  (5.2). 


dt 


An+0]an(O-^^y(O+^WO-«nW-«oW-yW) 


^M=-x0a0(tye0^r(t)+<iKxn(t) 


(5.5) 


(5.6) 


dy(t)  = 

dt 


enan(t)+e0a0(t) 


m 


y(0  +  Xnan  (0  +  A0a0  (f )  +  [S0a0  (t)  +  Snan  (t)][l  -  T)(t)]  (5.7) 


^  =  [{„«„(0  +  &«oW][l-*)] 


#M  _ 


dt 


=  5 


(5.8) 
(5.9) 


6.  Stochastic  Differential  Equation  Models 

The  size  of  the  state  space  makes  the  time-dependent  behavior  of  the 
continuous-time  Markov  chain  models  of  Sections  3-5  difficult  to  study.  One 
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approach  to  studying  the  behavior  is  by  Monte  Carlo  simulation;  cf.  Bois,  et  al. 
(1992).  Another  approach,  adopted  here,  is  to  approximate  the  continuous-time 
Markov  chain  model  by  a  diffusion  process.  Note  that  the  continuous-time 
Markov  chain  model  has  absorbing  states;  for  example,  any  state  with  Dn(t)  + 
D0(t)  =  0  is  absorbing.  Of  course  no  self-respecting  organ  would  ever  start  life  in 
such  a  state,  and  a  living  organ  would  presumably  die  from  other  unmodeled 
causes  long  before  such  a  state  is  reached,  i.e.  when  Dn(t)  +  D0(t)  <  D,  some  lower 
limit.  Barbour  (1976)  discusses  how  long  and  over  what  ranges  the  underlying 
continuous-time  Markov  chain  is  approximated  by  a  diffusion  process  of  the  type 
we  derive;  see  also  discussion  to  McNeil  and  Schach  (1973).  Stochastic 
differential  equation  models  can  be  written  for  all  the  models  of  Sections  3-5. 
We  will  illustrate  our  approach  by  writing  down  a  system  of  stochastic 
differential  equations  for  the  multivariate  process  (Dn(t),  D0(t),  G(t),  C(t),  T(0)  for 
Markov  generator  (5.2).  The  system  of  stochastic  differential  equations  is  as 
follows;  as  explained  subsequently,  the  dW-terms  are  normally  distributed. 


dDn(t)  =  -4(T(t))Dn(l)d(-e;(r(())G(()(D„(()/C(0>it 

-i(T(t))Dn(t)dt  +  M*(T(0)(C(0  -  D„(()  -  D0(()  -  G(())d(. 


drift = conditional 
expected  change 
m(t,t+dt) 


-JeQT^DMdW^  (t)  -  ^(T(t))Dn(t)dW^(t) 
-^ee:(T(t))G(t)(Dn(t)/C(t))dWen  (t) 

+^  (T(t))(C(t)  -  DH[t)  -  D0(t)  -  G(f))i  dW^t) 


f  dif  fusion = random  A 
change  in  (t,t  +  dt)  ) 


(6.1) 
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(drift) 


dD0(t)  =  -£(T(t))D0(t)dt  -  6'0(T(t))G(t)(D0(t)/C(t))dt 
+i(T(t))Dn(t)dt 

-JeX*0{T(t))D0(t)dWXo  (t)  +  ^(T(t))Dn(t)dW^(t) 
-pl{T(t))G(t){D0(t)/C(t))dW6o  (t) 


►(diffusion) 


(6.2) 


dG(t)  =  -e*n(T(t))G(t)(Dn(t)/C(t))dt  -  el(T(t))G(t)(D0(t)/C(t))dt 
+%(T(t))D0(t)dt+£tt(T(t))Dn(t)dt 


H*n(T{t))Dn(t) 


1- 


C(t) 


c, 


o . 


+ 


Simp,® 


i- 


Co. 


(drift) 


-^ee'n{T(t))G(t)(Dn(t)/C(t))dWen  (() 
-p'0(T(t))G(t)(D0(t)/C(t))dWec  (<) 
+^K(T('))Do(<)  ^  (')  +  ^(TOPnWW^  (() 


+   ^(r(())D„(() 


1- 


coJ 


dw|n(() 


+,K(T(0)Do(0 


C(() 


-0  J 


dw|o(f) 


(diffusion) 


(6.3) 
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mm&fm^-^s 


+&{T(t))Dn(t) 


C0J 


►(drift) 


+  e£(T(t)pB(t) 


1- 


co  J 


dW^t) 


+  e£(T(t)pn(t) 


r_£W 


«W£,(0 


(diffusion) 


(6.4) 


^  =  ^(T(0,Dn(0,Do(0,T*(0^}(drift) 

+aT{T(t)>Dn(t)/D0(t)fr(t))dWaT(t) 
+y[C^ardWT(t) 


(6.5) 


(diffusion) 


where   {w^^W},^ 

and  (w^f)}  are  independent  standard  Brownian  motions;  each  dW(0-term  is 
thus  Gaussian  with  mean  zero  and  variance  dt. 

Setting  the  constant  e  =  1  lets  the  variances  of  the  change  in  the  number  of 
active  cells,  ghosts  in  identified  states,  and  active  spaces  be  equal  to  the  mean 
respective  changes;  this  represents  Poisson  variability.  Setting  the  constant  £  >  1 
implies  that  the  variability  is  larger  than  Poisson  variability.  This  additional 
variability  may  be  the  result  of  inhomogeneities  in  the  organ  that  are  not 
explicitly  modeled;  cf.  Bass,  Robinson  and  Bracken  (1978)  for  discussion  of  a 
distributed  liver  model  with  random  inhomogeneities.  Such  over  dispersion  is 
frequently  encountered  in  practice,  and  may  be  accounted  for  parsimoniously  as 
we  do;  see  McCullagh  and  Nelder  (1983). 

We  assume  that  as  Co  — >  °° 


c0 


Co 


Co 
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in  probability  and  rewrite  Dn(t),  D0(t),  G(t),  C(t),  and  7X0  as  follows. 

Dn(t)  =  C0an(t)  +  Jq;xn(t)  (6.7a) 

D0(t)  =  Coa0(t)  +  ^X0(t)  (6.7b) 

'   G(0  =  0)7(0 +  VQ)xg(0  <6-7c> 

C^^Co^  +  ^X^t)  (6.7d) 

T(f)  =  Coi3(0  +  ^OyW  (67e) 

<***(*)  =  Cqt(0  +  ^ardWT(t)  (67f) 
The  transition  rates  will  be  scaled  and  expanded  as  follows: 

4(r(0)= A„(r(f)/c0)= a,0W))+ KiPCW')/^ +0(1)  (6.8a) 

A;(T(())  =  ^(T(0/Q,)  =  A0(/5(<))  +  Ji(«0)*(0A/5>  +  0(1)  (6.8b) 

/i*(r(t)) = m(t(0/co) = mWO) + nW))n)l-fo + o(i)  (6.8c) 

<£(T(0) = e„(T(t)/c0) = e„(^(<))+ ^(^(())v(()/VQ  +  o(i)  (6.8d) 

0'o{T(>))  =  *„C(0/C6)  =  *<,(/5(0)  +  %(m)n')/^  +  O(l)  (6.8e) 

«;(T(0)  =  ^(T(0/C0)  =  fe,(«0)+«;(/((0)V(0/>/Q+O(l)  (6.80 

&*(r(0) = lo(r(()/c0) = S,  WO) + £  (/»(0)*(0A/Q + o(i)  (6.8g) 
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$*(D„(t),Do(f),T(0,  %  (*))  =  c05(dm(0/c0  ,do(0/q,  ,t(0/c0  ,  t  (O/Co) 

=  Cb5(aB(f),ao(0,/5(*M0)  (6.8h) 

+Co5i  Siffi  +  Cb^i  ^  + 0^^  +  0(1) 

ar(DB(0,Do(0,r(0,T*(0)  =  VCo^r(D„W/C0  ,D0(f)/Q>  ,T(0/Q> ,  t  (0/q,) 

(6.8i) 

=  VQar(an(0,ao(0^(0^W)  +  0(1) 
Dividing  equations  (6.1)  -  (6.5)  by  Co,  scaling  and  substituting  (6.7a)  -  (6.7 f) 
results  in  a  system  of  equations  involving  the  deterministic  components  an(t), 
«o(0,  )0),  77(0,  and  p(t)  and  the  corresponding  stochastic  components  Xn(t),  X0(t), 
Xg(t),  Xc(0,  and  Y(t).  See  Appendix  1  for  details. 

6.1  Second  Moments 

The  stochastic  differential  equations  (A.l)  -  (A.10)  of  Appendix  1,  derived  for 
the  stochastic  disturbances,  can  be  written  in  matrix-vector  form  as 

dZ(t)  =  Aft)Z(f)  +  B(0  dmt)  (6.9) 

where  in  the  present  example  Z  is  a  5x1  column  vector,  A  is  a  5x5  matrix,  dW(t)  is 
10x1  column  of  independent  Gaussian  noise  terms  and  B(0  is  a  5x10  matrix. 
It  follows  from  theorem  (8.5.5)  of  Arnold  (1974)  that 

d 


|-E[z?(0]  =  2XAi;(0E[zi(0Z;(0]  +  l(B(/W)  (6.10) 


dt 


and 
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^E[z,(()z,(0]=lAft(f)4zt(0z,W] 

k 

+SA/Jt(OE[Z*(OZ,.(0]  (6.11) 

k 

k 
Thus,the  approximate  distribution  of  Dn(t),  (respectively  D0(t),  G(t),  C(t),T(t))  is 

normal  with  mean  Coan(t),  (respectively  Coctoit),  Coy(t)  CnrKO,  Cop(t))  and  variance 

C0e[x2(0],    (respectively    C0E[xo2(0],C0E[4(f)],CoE[xc2(0],CoE[v2(0])-    The 
approximate  covariance  of  Dn(t)  and  D0(t)  is  CoE[Xn(t)X0(t)],  etc. 

We  point  out  that  exactly  equivalent  results  can  be  obtained  by  introducing 
transforms  for  the  state  variables,  computed  by  taking  conditional  mathematical 
expectations  of  exponential  functions  of  the  Markov  generators.  This  is  the  route 
followed  by  McNeil  and  Schach  (1973);  see  also  Carpenter,  Gaver,  and  Jacobs 
(1993).  The  transform  approach  provides  a  verification  that  the  asymptotics 
described  here  actually  lead  to  limiting  Gaussian  distributions  and  to  Ornstein- 
Uhlenbeck  processes.  See  Ethier  and  Kurtz  (1986)  for  a  more  rigorous  discussion. 
The  approach  taken  here,  while  heuristic,  is  intuitively  appealing  and  reaches  the 
same  ultimate  conclusion. 

6.2  Numerical  Examples 

In  this  section  we  present  three  numerical  examples. 
a.   Constant  Rate  of  Input  of  Toxin 

We  assume  chronic  dosage  of  the  organ  by  toxin  that  is  delivered  there  at 
constant  mean  rate  but  with  substantial  randomness  around  that  mean. 

Figures  1-5  depict  the  time  variation  of  the  expected  numbers  of  new  cells, 
of  old  cells,  of  total  number  of  active  (new  plus  old)  cells,  number  of  signals 
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(ghosts),  and  amount  of  toxin  plus/minus  two  population  standard  deviations. 
There  is  no  organ  growth.  The  initial  conditions  are  the  limiting  moments  with 
no  toxin  input.  The  toxin  input  rate  0.2  units/cell  for  all  time.  The  new  cells 
without  toxin  have  mean  time  100  before  transitioning  to  old.  The  old  cells  have 
mean  lifetime  of  100  without  toxin.  The  number  of  cells  in  the  organ  is  1.5x10s. 
The  organ  size  is  about  that  of  a  mouse  liver.  Since  the  cells  in  a  mouse  liver  have 
a  mean  lifetime  of  about  200  days;  cf.  Bois,  et  al.  (1992),  the  above  is  consistent 
with  mouse  liver  behavior.  The  MATLAB  4th  and  5th  order  Runge-Kutta- 
Fehlberg  numerical  integration  method  was  used  to  obtain  the  solutions.  The 
MATLAB  plotting  algorithms  were  used  to  produce  the  figures. 

There  is  extra  Poisson  variability:  £  =  104. 

Note  that  the  total  mean  number  of  active  (new  and  old)  cells  decreases  as  the 
result  of  toxin  input,  as  seems  intuitively  reasonable.  The  mean  number  of  new 
cells  initially  decreases,  but  eventually  increases  to  a  new  steady-state  value 
which  is  larger  than  the  value  with  no  toxin.  The  mean  number  of  old  cells 
initially  increases  briefly,  but  then  decreases  to  a  new  steady-state  value  below 
the  value  for  no  toxin.  The  mean  amount  of  toxin  and  the  mean  number  of 
signals  (ghosts)  both  increase  to  new  steady-state  values. 

The  above  hypothetical  example  illustrates  how  model  parametric  inputs 
translate  into  an  account  of  the  dynamics  of  subpopulations  of  cells.  The  details 
of  the  transition  from  one  steady-state  situation  to  another  are  of  possible  interest 
in  that  the  organ  may  be  in  jeopardy  during  that  transitional  period  because  of, 
for  example,  new  cell  initial  downward  fluctuation.  The  effect  on  organ  mortality 
is,  however,  not  explicitly  modeled. 
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b.  A  Toxin  Pulse  or  Bolus  Dose 

Figures  6-10  display  the  mean  number  of  new  cells,  old  cells,  active  cells 
(new  and  old),  signals  (ghosts)  and  amount  of  toxin  plus  and  minus  2  standard 
deviations  for  the  case  in  which  the  rate  of  toxin  input  is  0.2  per  cell  for  the  first 
50  time  units  and  0  thereafter.  Once  again  there  is  no  organ  growth.  The  initial 
conditions  are  steady-state  values  when  there  is  no  toxin  input.  Note  that  the 
number  of  active  cells  decreases  during  toxin  input,  then  increases  and  slightly 
overshoots  the  steady-state  value  with  no  toxin  before  returning  to  its  steady- 
state  value  with  no  toxin. 

c.  Organ  Growth 

The  model  of  Section  5.1  and  (6.1)  -  (6.5)  is  initialized  with  the  steady-state 
values  for  the  second  moments  with  no  input  of  toxin.  However,  only  1/2  of  the 
available  spaces  are  active  (may  contain  cells);  the  steady-state  mean  of  new  and 
old  cells  and  ghosts  for  no  toxin  input  are  multiplied  by  1/2  and  used  as  initial 
values.  Figures  11-15  display  the  mean  numbers  of  new,  old,  and  active  cells, 
ghosts,  and  active  spaces  pi  us /minus  two  standard  deviations  for  the  case  of  no 
toxin  input  and  parameters  £n/o  =  0.5,  £n/i  =  0.5,  <f^o  =  1/  <5o,l  =  0.5.  The  mean 
number  of  active  spaces  is  the  maximum  number  by  time  20.  However,  the  mean 
number  of  old  and  new  cells  is  still  adjusting  by  time  100.  The  mean  number  of 
old  cells  initially  decreases  before  increasing.  The  mean  number  of  new  cells 
initially  increases  before  decreasing. 

The  present  model  might  be  useful  for  describing  the  aftermath  of  a  partial 
hepatechtomy:  about  one-half  of  a  mature  liver  remains,  and  grows  back  with 
subpopulations  of  cells  responding  as  shown.  Note  that,  according  to  the  present 
model  in  Figure  12,  the  organ  must  withstand  an  early  signal  to  old  cells  to 
undergo  substantial  depletion  to  enter  mitosis  and  bring  forth  replacements;  the 


27 


latter  are  the  new  cells  depicted  in  Figure  11.  It  may  be  that  the  sudden  depletion 
of  old  cells  to  their  low  point  at  about  t  -  6  (Figure  12)  actually  puts  the  organ  in 
jeopardy. 

Figures  16  -  21  display  the  mean  number  of  new,  old,  and  active  cells,  ghosts 
and  active  spaces  plus  /minus  two  standard  deviations  for  the  same  parameters 
as  used  in  Figures  11-15  but  with  two  levels  of  toxin  input.  The  dashed  lines 
correspond  to  no  toxin  input.  The  solid  lines  correspond  to  a  constant  toxin  input 
of  0.2  per  cell  per  unit  time.  The  graphs  with  positive  toxin  input  have  the  same 
general  shape  as  those  with  no  toxin  input.  The  mean  number  of  new  cells  is 
initially  less  for  a  positive  toxin  input  than  for  no  toxin  input.  However,  by  time 
100  the  mean  number  of  new  cells  with  positive  toxin  input  is  becoming  larger 
than  that  for  no  toxin  input. 

The  effect  of  positive  toxin  input  on  the  mean  number  of  old  cells  is  that  it  is 
initially  larger  than  that  for  no  toxin  but  by  time  100  the  mean  number  of  old 
cells  with  positive  toxin  is  less  than  that  for  no  toxin.  The  effect  of  positive  toxin 
input  on  the  mean  number  of  active  (old  and  new)  cells  is  to  decrease  it  and 
delay  its  approach  to  a  steady  state  value.  The  mean  number  of  ghosts  is  larger 
with  positive  toxin  input.  Positive  toxin  input  delays  the  overall  growth  in  the 
number  of  active  spaces. 
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7     0.5  of  the  spaces  are  active 
x  10       mean  number  or  old  cells  plus  minus  2  standard  deviations 
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0.5  of  the  spaces  are  active  initially 

mean  number  of  ghosts  plus  minus  2  standard  deviations 
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Appendix  1 

This  appendix  provides  the  detailed  asymptotic  development  of  both  the 
deterministic  approximation  (5.5)  -  (5.9)  but  also  for  the  stochastic  components 
introduced  in  (6.7).  Start  by  introducing  the  normalization  (6.7)  into  (6.1)  -  (6.5) 
and  divide  by  Co-  The  following  results  are  obtained: 
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The  terms  of  order  1  in  (A.l)  -  (A.5)  yield  the  deterministic  equations  (5.5) 
(5.9). 
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The  terms  of  order  1/VCo  in  (A.l)  -  (A.5)  give  the  following  stochastic 
differential  equations  for  the  noise  terms. 


-e„(«0)x,(()^-^(«/))rW^v(() 


-»»(W)M0 


X„(()-^Xc(f) 


m 


w 


-<p{P(t))Xn(t)-<p'(p(t))an(t)Y(t) 

MP(t)pc(t)-x„(t)-x0(t)-xg(t)] 
+M'(/3(0)[*)-«n(')-«o(0-r(0M<) 


-^„(«t))r(0^dw,n(() 

+^{P(t)){v(t)  -an(t)-a0(t)-  r(())  d  W„  (() 


(A.6) 


54 


BM = -a0(/j(o)x0(()  -  ^(wKcmo 


-».(«t)MO 


-Lxo(()-^xc(0 
w)         i(0 


+0(/}(<))X„(f)  +  0'(/5(')K(W<) 


-^(j8(i)KW^  (f ) + ^(/5(')KWdvv0(() 


(A.7) 


-U(/J(0M0^!<*wflott 


55 


^  =  -eMt))xg(t)^-emm^y(t) 


-on(P(t))y(t) 


xn(t)-^xc(t) 


m 


m 


Mm)x8(t)^-e'Mt))Y(t)^Y(t) 


-e0{P(t))r(t) 


J-X0(t)-^Xc(t) 
W)  77(f) 


-zn(P(t)Kmc(t) 

Hn  («0)(i  -  n(t))xn  (t) + s (y3(0K(0(i  -  u(0)*W 
-t(W))ao(0xc(0+6,(«0Xi-*))XoW 

+^(/?(0)ao(0(l-r?(0M0 
+Ao(^(0)Xo(0  +  A^(«0K(0V(0 

-V^n  («OM'K(OA?M  <*w*n  (0  -  M(«0)r  (0«o  (OMO  ^\  (0 
+<J&oWt))c*o(t)  dwXo  (0 + ^(«t)KW%  (0 

+^„(«0K  (t)[l  - 1,(0]  ^  (0  +  JeS0(P(t))<*o(t  jp  "  *)]  Wfi,  (0 


(A.8) 


56 


&j&  =  -£,  (P(t))a0(t)Xc(t)  +  sW»K«i-*)W*) 

+«o(/?(0)[l-77(0]Xo(0 
+{„(«0)[l-r?(0]Xn(0 

+Vfifc(«OKW[1-»KO]^  (0 

^  =  ^(^(f)/an(f),ao(0,rW)Y(0 
+5in(j3W/an(0/Ofo(f)/T(0)Xn(0 

+cT(p(t)^n(t)foc0(t)f<t))dWaT(t) 


(A.9) 


(A.10) 


+crTrfWT(f) 


57 


REFERENCES 

L.  Arnold,  Stochastic  Differential  Equations:  Theory  and  Applications.  J.  Wiley  and 
Sons,  New  York,  1974. 

A.  D.  Barbour,  "Quasi-stationary  distributions  in  Markov  population  processes," 
Adv.  Appl.  Prob.,  8  (1976),  pp.  296-314. 

L.  Bass,  P.  Robinson,  and  A.  J.  Bracken,  "Hepatic  elimination  of  flowing 
substrates:  the  distributed  model,"  /.  Theoretical  Biology,  72  (1978),  pp.  161- 
184. 

F.  Y.  Bois  and  P.  J.  E.  Compton-Quintana,  "Sensitivity  analysis  of  a  new  model  of 
carcinogenesis."  /.  Theor.  Biol,  159  (1992),  pp.  361-375. 

R.  L.  Carpenter,  D.  P.  Gaver,  and  P.  A.  Jacobs,  "An  exploratory  stochastic  model 
for  toxic  effects  on  cells,"  Naval  Postgraduate  School  Technical  Report,  NPS- 
OR-93-014,  (1993),  pp.  1-76. 

S.  N.  Ethier  and  T.  G.  Kurtz,  Markov  Processes:  Characterization  and  Convergence, 
John  Wiley  and  Sons,  1986. 

H.  I.  Freedman  and  J.  B.  Shukla,  "Models  for  the  effect  of  toxicant  in  single- 
species  and  predator-prey  systems."  /.  Math  Biology,  30  (1991),  pp.  15-30. 

P.  Jagers,  "Stochastic  models  for  cell  kinetics."  Bulletin  of  Math.  Biology,  45,  No.  4 
(1983),  pp.  507-519. 

R.  L.  Kodell,  D.  Krewski,  and  J.  M.  Zielinski,  "Additive  and  multiplicative 
relative  risk  in  the  two-stage  clonal  expansion  model  of  carcinogenesis," 
Risk  Analysis,  11  (1991),  pp.  483-490. 

D.  Krewski  and  C.  Franklin  (eds.)  Statistics  in  Toxicology,  Gordon  and  Breach 
Science  Publishers,  New  York,  1991. 

J.  P.  Lehoczky  and  D.  P.  Gaver,  "A  diffusion-approximation  analysis  of  a  general 
n-compartment  system,"  Mathematical  Biosciences,  36  (1977),  pp.  127-148. 


58 


The  Math  Works,  Inc.  MATLAB  Reference  Guide,  The  Math  Works,  Inc.,  Natick, 
MA,  August  1992. 

P.  McCullagh  and  J.  A.  Nelder,  Generalized  Linear  Models,  Chapman  and  Hall  Ltd., 
New  York,  1983. 

D.  R.  McNeil  and  S.  Schach,  "Central  limit  analogues  for  Markov  population 
processes,"  (with  discussion),  Journal  of  the  Royal  Statistical  Soc,  1  (1973), 
pp.  1-23. 

S.  H.  Moolgavkar,  A.  Dewanji,  and  D.  J.  Venzon,  "A  stochastic  two-stage  model 
for  cancer  risk  assessment,  I:  the  hazard  function  and  the  probability  of 
tumor."  Risk  Analysis,  8  (1988),  pp.  383-392. 

P.  J.  Robinson,  A.  N.  Pettitt,  J.  Zornig,  and  L.  Bass,  "A  Bayesian  analysis  of 
capillary  heterogeneity  in  the  intact  pig  liver,"  Biometrics,  39  (1983),  pp.  61- 
69. 

L.  A.  Segel  and  M.  Slemrod,  'The  quasi-steady-state  assumption:  a  case  study  in 
perturbation."  SI  AM  Review,  31,  No.  3,  (1989),  pp.  446-477. 


59 


INITIAL  DISTRIBUTION  LIST 


1.  Research  Office  (Code  08) 1 

Naval  Postgraduate  School 

Monterey,  C A   93943-5000 

2.  Dudley  Knox  Library  (Code  52) 2 

Naval  Postgraduate  School 

Monterey,  CA   93943-5002 

3.  Defense  Technical  Information  Center 2 

Cameron  Station 

Alexandria,  VA   22314 

4.  Department  of  Operations  Research  (Code  OR/Bi) 1 

Naval  Postgraduate  School 

Monterey,  CA   93943-5000 

5.  Prof.  Donald  P.  Gaver  (Code  OR/Gv) 10 

Naval  Postgraduate  School 

Monterey,  CA   93943-5000 

6.  Prof.  Patricia  A.  Jacobs  (Code  OR/Jc) 5 

Naval  Postgraduate  School 

Monterey,  CA   93943-5000 

7.  Dr.  J.  Abrahams,  Code  1111,  Room  607 1 

Mathematical  Sciences  Division,  Office  of  Naval  Research 

800  North  Quincy  Street 
Arlington,  VA   22217-5000 

8.  Dr.  John  Bailar 1 

468NSt.NW 
Washington,  DC  20024 

9.  Dr.  David  Brillinger 1 

Statistics  Department 

University  of  California 
Berkeley,  CA   94720 


60 


10.     Dr.  David  Burman 

AT&T  Bell  Telephone  Laboratories 
600  Mountain  Avenue 
Murray  Hill,  NJ    07974 


11.    Prof.  Brad  Carlin 

School  of  Public  Health 
University  of  Minnesota 
Mayo  Bldg.  A460 
Minneapolis,  MN   55455 


12.  Dr.  Robert  Carpenter 

NAMRI/Navy  Toxicology  Detachment 

Bldg.  433,  Area  B 

Wright-Patterson  AFB,  OH  45433-6503 

13.  Center  for  Naval  Analyses 

4401  Ford  Avenue 

Alexandria,  VA   22302-0268 

14.  Dr.  Edward  G.  Coffman,  Jr 

AT&T  Bell  Telephone  Laboratories 
600  Mountain  Avenue 

Murray  Hill,  NJ   07974 

15.  Professor  Sir  David  Cox 

Nuffield  College 

Oxford,  OXI  INF 
ENGLAND 


16.    Professor  H.  G.  Daellenbach 

Department  of  Operations  Research 
University  of  Canterbury 
Christchurch,  NEW  ZEALAND 


17.    Dr.  D.  F.  Daley 

Statistic  Dept.  (I.A.S.) 
Australian  National  University 
Canberra,  A.C.T.   2606 
AUSTRALIA 


18.    Dr.  Naihua  Duan 

RAND  Corporation 
Santa  Monica,  CA  90406 


61 


19.    Prof.  Bradley  Efron... 
Statistics  Dept. 
Sequoia  Hall 
Stanford  University 
Stanford,  CA   94305 


20.    Prof.  George  S.  Fishman 

Curr.  in  OR  &  Systems  Analysis 
University  of  North  Carolina 
Chapel  Hill,  NC   20742 


21.    Henry  S.  Gardner 

U.S.  Army  Biological  Research  &  Development  Laboratory 
Ft.  Detrick,  MD  21702-5010 


22.  Dr.  Andrew  Gelman 

Statistics  Dept. 
University  of  California 
Berkeley,  CA  94720 

23.  Dr.NeilGerr 

Office  of  Naval  Research 
Arlington,  V A   22217 


24.  Prof.  Peter  Glynn 

Dept.  of  Operations  Research 
Stanford  University 
Stanford,  CA   94350 

25.  Prof.  Linda  V.  Green 

Graduate  School  of  Business 
Columbia  University 

New  York,  NY   10027 

26.  Prof.  J.  Michael  Harrison 

Graduate  School  of  Business 
Stanford  University 
Stanford,  CA  94305-5015 


27.    Dr.  P.  Heidelberger 

IBM  Research  Laboratory 
Yorktown  Heights 
New  York,  NY   10598 


62 


28.  Dr.  D.  C.  Hoaglin 

Department  of  Statistics 
Harvard  University 

1  Oxford  Street 
Cambridge,  MA   02138 

29.  Prof.  D.  L.  Iglehart 

Dept.  of  Operations  Research 
Stanford  University 
Stanford,  CA   94350 

30.  Institute  for  Defense  Analysis 
1800  North  Beauregard 
Alexandria,  VA   22311 

31.  Prof.  J.  B.  Kadane 

Dept.  of  Statistics 
Carnegie-Mellon  University 
Pittsburgh,  PA    15213 


32.  Dr.  F.  P.  Kelly 

Statistics  Laboratory 

16  Mill  Lane 

Cambridge 

ENGLAND 

33.  Dr.  Jon  Kettenring 

Bellcore 

445  South  Street 

Morris  Township,  NJ   07962-1910 

34.  Koh  Peng  Kong 

OA  Branch,  DSO 

Ministry  of  Defense 
Blk  29  Middlesex  Road 
SINGAPORE    1024 


35.    Prof.  Guy  Latouche 

University  Libre  Bruxelles 
C.P.  212,  Blvd.  De  Triomphe 
B-1050  Bruxelles 
BELGIUM 


63 


36.    Dr.  A.  J.  Lawrence 

Dept.  of  Mathematics, 
University  of  Birmingham 
P.  O.  Box  363 
Birmingham  B15  2TT 
ENGLAND 


37.  Prof.  M.  Leadbetter 

Department  of  Statistics 
University  of  North  Carolina 
Chapel  Hill,  NC   27514 

38.  Prof.  J.  Lehoczky 

Department  of  Statistics 
Carnegie-Mellon  University 
Pittsburgh,  PA    15213 

39.  Dr.  Colin  Mallows 

AT&T  Bell  Telephone  Laboratories 
600  Mountain  Avenue 

Murray  Hill,NJ   07974 


40.    Dr.  Sati  Mazumdar 

Biostatistics  Dept. 
University  of  Pittsburgh 
Graduate  School  of  Public  Health 
Pittsburgh,  PA  15261 


41.  Dr.  James  McKenna 

Bell  Communications  Research 
445  South  Street 
Morristown,  NJ   07960-1910 

42.  Prof.  Carl  N.  Morris 

Statistics  Department 
Harvard  University 

1  Oxford  St. 
Cambridge,  MA   02138 


43.    Dr.  John  A.  Morrison 

AT&T  Bell  Telephone  Laboratories 
600  Mountain  Avenue 
Murray  Hill,  NJ   07974 


64 


44.  Prof.  F.  W.  Mosteller 

Department  of  Statistics 

Harvard  University 
1  Oxford  St. 
Cambridge,  MA   02138 

45.  Operations  Research  Center,  Rm  E40-164 
Massachusetts  Institute  of  Technology 
Attn:  R.  C.  Larson  and  J.  F.  Shapiro 
Cambridge,  MA   02139 


46.     Dr.  John  Orav 

Biostatistics  Department 
Harvard  School  of  Public  Health 
677  Huntington  Ave. 
Boston,  MA    02115 


47.    Dr.  Jim  Petty 

National  Biological  Survey 
4200  New  Haven  Road 
Columbia,  MO  65201 


48.  Dr.  V.  Ramaswami 

MRE  2Q-358 

Bell  Communications  Research,  Inc. 
445  South  Street 
Morristown,  NJ   07960 

49.  Dr.  Martin  Reiman 

Rm#2C-117 

AT&T  Bell  labs 
600  Mountain  Ave. 
Murray  Hill,  NJ    07974-2040 


50.  Dr.  Lorenz  Rhomberg 

Harvard  Center  for  Risk  Analysis 
Harvard  University 
Cambridge,  MA  02138 

51.  Prof.  Maria  Rieders 

Dept.  of  Industrial  Eng. 
Northwestern  Univ. 
Evanston,IL   60208 


65 


52.    Dr.  John  E.  Rolph 

RAND  Corporation 

1700  Main  St. 

Santa  Monica,  CA   90406 


53.    Prof.  Frank  Samaniego  .. 
Statistics  Department 
University  of  California 
Davis,  CA   95616 


54.    Prof.  W.  R.  Schucany 

Dept.  of  Statistics 

Southern  Methodist  University 

Dallas,  TX   75222 


55.    Prof.  D.  C.  Siegmund 
Dept.  of  Statistics 
Sequoia  Hall 
Stanford  University 
Stanford,  CA   94305 


56.  Prof.  H.  Solomon 

Department  of  Statistics 

Sequoia  Hall 
Stanford  University 
Stanford,  CA   94305 

57.  Andrew  Solow 

Woods  Hole  Oceanographic  Institute 
Woods  Hole,  MA  02543 


58.  Prof.  W.  Stuetzle 

Department  of  Statistics 
University  of  Washington 

Seattle,  WA   98195 

59.  Dr.  D.  J.  Svensgaard 

Health  Effects  Research  Lab 

U.  S.  Environmental  Protection  Agency 
Biostatistics  Branch  MD55 
Research  Triangle  Park,  NC  27711 


60.    Prof.  J.  R.  Thompson 

Dept.  of  Mathematical  Science 
Rice  University 
Houston,  TX   77001 


66 


61.  Prof.  J.  W.  Tukey 

Statistics  Dept.,  Fine  Hall 

Princeton  University 
Princeton,  NJ   08540 

62.  Dr.  D.  Vere-Jones 

Dept.  of  Math,  Victoria  Univ.  of  Wellington 
P.  O.  Box  196 

Wellington 
NEW  ZEALAND 


63.    Prof.  David  L.  Wallace 

Statistics  Dept.,  University  of  Chicago 
5734  S.  University  Ave. 
Chicago,  IL   60637 


64.  Dr.  Ed  Wegman 

George  Mason  University 
Fairfax,  VA  22030 

65.  Dr.  Alan  Weiss 

Rm.  2C-118 

AT&T  Bell  Laboratories 
600  Mountain  Avenue 
Murray  Hill,  NJ   07974-2040 


66.  Dr.  P.  Welch 

IBM  Research  Laboratory 
Yorktown  Heights,  NY   10598 

67.  Prof.  Roy  Welsch 

Sloan  School 

M.I.T. 

Cambridge,  MA   02139 


68.     Dr.  Tandore  Narayanan 

NAMRI/Navy  Toxicology  Detachment 

Bldg.  433,  Area  B 

Wright-Patterson  AFB,  OH  45433-6503 


67 


DUDLEY  KNOX  LIBRARY 


3  2768  00331 674  6 


